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ABSTRACT 

Upcoming weak lensing surveys, such as LSST, EUCLID, and WFIRST, aim to mea- 
sure the matter power spectrum with unprecedented accuracy. In order to fully exploit 
these observations, models are needed that, given a set of cosmological parameters, 
can predict the non-linear matter power spectrum at the level of 1% or better for 
scales corresponding to comoving wave numbers 0.1 ^ fc ^ 10/iMpc~^. We have em- 
ployed the large suite of simulations from the OWLS project to investigate the effects 
of various baryonic processes on the matter power spectrum. In addition, we have 
examined the distribution of power over different mass components, the back-reaction 
of the baryons on the CDM, and the evolution of the dominant effects on the matter 
power spectrum. We find that single baryonic processes are capable of changing the 
power spectrum by up to several tens of per cent. Our simulation that includes AGN 
feedback, which we consider to be our most realistic simulation as, unlike those used in 
previous studies, it has been shown to solve the overcooling problem and to reproduce 
optical and X-ray observations of groups of galaxies, predicts a decrease in power rel- 
ative to a dark matter only simulation ranging, at z = 0, from 1% at fc f» 0.3 /iMpc~^ 
to 10% at A; « 1 /iMpc^^ and to 30% at fc w 10/iMpc~^. This contradicts the naive 
view that baryons raise the power through cooling, which is the dominant effect only 
for k > 70/iMpc~^. Therefore, baryons, and particularly AGN feedback, cannot be 
ignored in theoretical power spectra for k > 0.3/iMpc~^. It will thus be necessary 
to improve our understanding of feedback processes in galaxy formation, or at least 
to constrain them through auxiliary observations, before we can fulfil the goals of 
upcoming weak lensing surveys. 

Key words: Cosmology: theory, large-scale structure of Universe, galaxies: formation, 
gravitational lensing: weak, surveys 



1 INTRODUCTION 

One of the aims of cosmology is to find the initial conditions 
for structure formation in the Universe. These can be char- 
acterised by a single set of cosmological parameters, which 
directly influence the formation, growth and clustering of 
structure, and hence the distribution of matter as we ob- 
serve it today. 

A powerful measure of the statistical distribution of 
matter (and a sufficient one for the case of Gaussian fluctu- 
ations), is the matter power spectrum, P{k), where k is the 
comoving wave number corresponding to a comoving spatial 
scale A — 2n/k. Given a sufficiently accurate model for the 
formation of structure, we can infer the initial, linear power 
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spectrum from the observed, non-linear one. Moreover, as 
the rate of growth of structure depends on the expansion 
history, such a model also allows us to convert observations 
of the evolution of the power spectrum into measurements of 
other cosmological parameters such as the equation of state 
of the dark energy. 

Some of the most accurate measurements of the mat- 
ter power spectrum come from studies of weak, gravi- 



tational lensing (e.g. iMassev et all l2007l : IFu et all I2OO8I: 



^ Schrabback et al]|2010l ). galaxy clustering (e.g. I Cole et al 



Reid et al.ll2010l) and the Lya forest (e.g. IViel et al 



McDonald et al.ll2006l ). Up to a few years ago, the 
statistical errors were sufficiently large that one could 
use analytical predictions (always assuming, amongst other 
things, that t he Universe contai n s only d ark matte r) , such 
as those by iPeacock fc Doddsl (|l996l ). iMa et all (|l999l ) 
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and ISmith et al.l ll2q03h. The latter used ideas from the 
"halo model" (e.g. IPeacock fc SmithI I2OO0I : ISeliakI I2OO0I : 
ICoorav fc Shethll2002l ) to improve upon the accuracy of sim- 
pler analytical predictions. In recent years the further im- 
provement of this model has become increasingly depen- 
dent on the results from A'^-body simulations, such as the 
deri ved concentration-mass relation for dark matter haloe s 
(e.g. iNeto et al.l[2007l : iDuffv et al.l[2008l : iHilbert et aLll2009h . 
If baryonic effects were negligible, then these methods would 
allow the matter power spectrum to be predicted with 
an accuracy of ~ 1% for wave numbers k < 1/iMpc^^ 
iHeitmann et al.ll201Cf ). However, we will show here that 
baryonic effects are larger than this on the scales relevant 
for many observations. 

Upcoming weak lensing surveys aim to measure the 
matter power spectrum on scales of 0.1/iMpc~^ < k < 
10/iMpc~^. In order to reach the level of precision their 
instruments are capable of, surveys such as LSStQ EU- 
CLID and WFIRSlH need to be calibrated using the- 
oretical models that retain 1% a,ccura cies on these scales 
IIHuterer fc Takadallioosi : iLaureiidlioO^ l This is, however, 
not as straightforward as increasing the resolution of exist- 
ing A'^-body simulations: many authors have demonstrated 
that on these scales baryonic matter, which is not ac- 
counted for in currently employ ed theoretica l models, intro- 
duce s deviations of up to 10% (lwhitd|2004l: IZhan fc Kii^ 
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Casarini et al. 2011^ 



2008 



Guillet et al.ll2010l : 



Recent hydrodynamic simulations include many of the 
physical processes associated with baryons, such as radiative 
cooling, star formation and supernova (SN) feedback. How- 
ever, the processes which cannot be resolved in simulations 
are generally also not entirely understood, and different pre- 
scriptions exist that aim to model the same physics. Because 
of this, different authors may find significantly different re- 
sults even when including the same baryonic processes. Fur- 
thermore, it is not a priori clear which physical effects are 
capable of changing the matter power spectrum at the 1% 
level and should therefore be included. These modelling un- 
certainties may thus prevent upcoming surveys from further 
constraining the cosmological parameters of our Universe. 

Here we employ a large suite of state-of-the-art cosmo- 
logical, hydrodyna mical simulations from the OWLS project 
l|Schave et all 12010 ) to systematically study the effects of 
various baryonic processes on the matter power spectrum 
over a wide range of scales, fc ~ 0.1 - 500/iMpc"^ These 
processes include metal-line cooling, different prescriptions 
for SN feedback, and feedback from active galactic nuclei 
(AGN). We will see that all of our results are heavily in- 
fluenced by the inclusion of AGN feedback, which was not 



^ http://www.lsst.org/lsst 

^ http://www.euclid-imaging.net/ 

^ http://wfirst.gsfc.nasa.gov/ 

* Since cosmological parameters are inferred from cosmic shear 
using a complicated weighting of the power spectrum over a range 
of scales and redshifts, the relation between the accuracy with 
which these parameters can be determined and the uncertainty 
in the models depends on the survey and is different for different 
parameters. Semboloni et al. (in preparation) will present a more 
detailed study of the consequences of our findings for weak lensing 
surveys. 



considered by earlier studies and which has been shown to 
solve the overcooling problem that has long plagued hydro- 
dynamical simulations and to lead to an excellent match to 
bot h the optical and X-ray pro perties of groups of galax;- 
ies ijMcCarthv et al.ll2010l . [20111 ). Outflows driven by AGN 
strongly increase the scale out to which baryons modify the 
power spectrum. We also investigate how the power is dis- 
tributed over different components (i.e. CDM, gas and stars) 
and examine the back-reaction of the baryons on the dark 
matter. In a follow-up paper (Semboloni et al., in prepa- 
ration), we will quantify the implications for current and 
proposed weak lensing surveys and we will show how the 
uncertainty due to baryonic physics can be reduced by mak- 
ing use of additional observations of groups and clusters. 

This paper is organised as follows. In i|2]we discuss the 
simulations and the power spectrum estimator employed. 
In our main results section, [JS] we compare our dark mat- 
ter only simulation to analytical estimates (SjSTT}, we com- 
pare power spectra of simulations with different baryonic 
processes ( ^3.2p . and we investigate how the power is dis- 
tributed over different physical components ( ii3.3p . In this 
section we also examine the back-reaction of galaxy forma- 
tion on the dark matter f ^3.4p and we consider the evolution 
of the most dominant effects on the power spectrum (|33]). 
We compare to the results found by other authors in SUJa^nd 
provide a summary in ij5] Finally, we test the convergence of 
our results in Appendix A and provide tables of the power 
spectra of all simulations in Appendix B. 

We note that all distances quoted in this paper are co- 
moving and all power spectra are obtained at redshift zero, 
unless stated otherwise. 



2 SIMULATIONS 



The OWLS project (|Schave et al.ll2010l ). where OWLS is 
an acronym for Overwhelmingly Large Simulations, is a 
suite of large, cosmological, hydrodynamical simulations. 
The code used is a heavily extended version of GADGET ill, a 
Lagrangian code which was last described in ISpringell()2005l l. 
It uses a TreePM algorithm to efficiently calculate the grav- 
itational forces (where PM stands for Particle Mesh and the 
"Tree" describes the structure in which the p articles are or- 
ganised fcjr this calculation , see for example iBarnes fc HutI 
|l986l : IXulll995l : lBaglall2002l ) and Smoothed Particle Hydro- 
dynamics (SPH) to follow and evolve the gas particles (see 
[Posswo g 2009 for a review). 

There are two main sets of simulations, which have pe- 
riodic boxes of size L = 25 and 100 comoving Mpc on 
a side, and are run down to redshifts z = 2 and 0, respec- 
tively. Most simulations use 512^ coUisionless cold dark mat- 
ter (CDM) particles and an equal number of baryonic (col- 
lisional gas or coUisionless star) particles. We will refer to 
the particle number used in a simulation with the parame- 
ter N = Afp^ft (= 512 for the high-resolution simulations). 
In this work we will focus on 2 = and hence on the sim- 
ulations using a 100 Mpc box. The particle masses are 
4.06 X 10® h-^ Mq [L/(100 h'^ Mpc)]^[Ar/512]-^ for the dark 
matter and 8.66 x lO'' h'^ M0[L/(1OO h'^ Mpc)f[N/^12]-^ 
for the baryons. The gravitational forces are softened on 
a comoving scale of 1/25 of the initial mean interparticle 
spacing, L/N, but the softening length is limited to a max- 
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imum physical scale oi 2h~^ kpc[I//(100 Mpc)] which is 
reached at z — 2.91. The SPH calculations use 48 neigh- 
bours. 

For the initial conditions, a theoretical matter power 
spectrum - which of course depends on the chosen set 
of cosmological parameters - is generated using cmbfast 
jSeliak fc Zaldarriagal [l99^ . version 4.1). Prior to imposing 
the linear input spectrum, the particle s are set up in an 
initially glass-like state, as described in IWhitd ()l994l ). The 
particles a re then evolved to redshift z — 127 using the 
IZel'Dovichl (|l97Q ) approximation. 

On small scales, the physics of gala:xy formation is un- 
resolved, and subgrid models are needed to include baryonic 
effects like radiative cooling, star formation and supernova 
feedback. Although each OWLS run is a state-of-the-art cos- 
mological simulation in itself, the real power of the OWLS 
project lies in the fact that it is composed of more than 
50 simulations that all incorporate different sets of physical 
processes, parameter values, or subgrid recipes. In this way 
the effects of turning off or tweaking a single process can be 
studied in detail, making it especially well-suited to investi- 
gate which processes can, by themselves, change the power 
at A; ~ 1 - 10/iMpc"^ by > 1%. In this paper we briefly 
describe the subgrid physics included in the reference simu- 
lation, as well as the differences with respect to simulations 
we compare to in 8^3.21 For a more detailed treatment of the 
simulat ions and the differe nt physics models included, we 
refer to lSchave et"al] (|2010l ). 



2.1 The reference simulation 

As the intention of the OWLS project is to investigate the 
effects of altering or adding a single physical process, it is 
convenient to have a single simulation that acts as the ba- 
sis for all other simulations. Such a "default" simulation 
should of course include many of the physical processes that 
we know to be important already, as ideally we would only 
want to vary one thing at a time. We call this simulation 
the reference simulation, or REF for short. Note that this 
is not intended to be the "best" simulation, but simply a 
model to build on. In fact, it has for example been shown 
that AGN feedback, which was not included in the REF 
model and which we briefly discuss in the next section, is 
required to match observations of g roups and clusters of 
galaxies ("M cCarthv et al.ll201(]| . 1201]] ). 

We assume cosmological parameter values derived from 
the Wilkins on Microwave Aniso tropy Probe (WMAP) 3- 
year results (|Spergel et al.llioOTi '): {Sim, f^b, SIa, erg, Ws, h} 
= {0.238, 0.0418, 0.762, 0.74, 0.951, 0.73}. Except for erg, 
all of these are cons istent with the WMAP 7-year data 
l|Komatsu et al.|[201ll ). This specific parameter describes the 
root mean square fluctuation in spheres with a radius of 
8h~^ Mpc linearly extrapolated to 2: = and effectively nor- 
malises the matter power spectrum. Measurements in the 
last few years have systematically increased the value of as, 
which may influence the validity of our results. To check the 
effects of using "wrong" values for this and other cosmo- 
logical parameters, we have re-run our two most important 
simulations - one with only dark matter and one in which 
AGN feedback is added to the reference model - using the 
WMAP7 cosmology. We briefly discuss these at the end of 



section i]2.2l As we shall see in ^3.51 this change in cosmology 
does not affect our conclusions. 

The reference simulation includes both radiative cool- 
ing and heating , which are modelled using the prescription of 
|w iersma et al.l (120091 '). Net radiative cooling rates are com- 
puted on an element-by-element basis i n the presence of 
the co smic microwave background and the lHaardt fc Madaul 
(|200lh model for the UV and X-ray background radiation 
from quasars and galaxies, taking into account the con- 
tributions of eleven different elements pre-computed using 
the publicly availab le photo-ioniz at ion p ackage CLOUDY, 
last described by iFerland et "all l| 19981 ). The effects of 
hydrogen ionizat i on ar e modelled by switching on the 
iHaardt fc Madaul l|200lh model at z = 9. 

Cosmological simulations do not yet come close to re- 
solving the process of star formation, and so a subgrid 
recipe has to be included for this as well. In our simu- 
lations, gas particles can be converted into star particles 
once their hydrogen number densities exceed the th reshold 
for th ermo-gravitational instability (tih = 0.1 cm~'^: ISchavj 
l2004h . Cold gas particles with higher densities follow an im- 
posed e2u^tion_of_statej_P« Here -ycS = 4/3, for 
which [Schave fc Dalla Vecchial l|2008h showed that both the 
Jeans mass and the ratio of the Jeans length to the SPH ker- 
nel are independent of the density, thus preventing spurious 
fragmentation due to a lack of numerical resolution. Using 
their pressure-dependent prescription for star formation, the 
observed Kennicutt-Schmidt relation, a surface density scal- 
ing law fo r the star forma tion rate that can be written as 
E« oc Eg (|Kennicuttl 1 19981 ) . is reproduced by construction, 
independent of the imposed equation of st ate. 

The reference simulation assumes a IChabrieil (|2003D 
stellar Initial Mass Function (IMF) with low and high mass 
cut-offs at 0.1 and 100 Mq, respectively. The release of 
hydrogen, helium and heavier elements by these stars to 
the surrounding gas is tracked as well: gas can be ejected 
through Type II SNe and stellar winds for massive stars, 
and Type la SNe and Asymptotic Giant Branch (AGB) 
stars for intermediate mass stars. This implementation of 
stellar evolution and chemical enrichment is discussed in 
IWiersma et"all (I2009D . 

Finally, the reference simulatio n includes a prescription 
fo r su pernova feedback, discussed in lDalla Vecchia fc Schavg 
( 2008) . Supernovae are capable of depositing a significant 
amount of energy in the surrounding gas, driving large- 
scale winds that may eject large amounts of gas, dramat- 
ically suppressing the formation of stars. In the model used 
here, the energy from SNe is injected into the gas kineti- 
cally. After a delay time of 30Myr, a new star particle j 
will "kick" a neighbouring SPH particle i with a probabil- 
ity r^rrij/ X]i=f'^ ^ random direction, giving it an extra 
velocity Uw The reference simulation uses the values rj = 2 
for the initial wind mass loading and Uw ~ 600 km s^^ for 
the initial wind velocity, which corresponds to 40% of the 
available kinetic energy for our IMF. 



2.2 Other models 

The OWLS project includes many variations on the refer- 
ence simulation. We will now briefly discuss the simulations 
that we compare to in H3.2I The different models are listed 
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Table 1. The different variations on the reference simulation that are compared in this paper. Unless noted otherwise, all simulations 
use a set of cosmological parameters derived from the WMAP3 results and use identical initial conditions. 



Simulation 


Description 


AGN 


Includes AGN (in addition to SN feedback) 


AGN.WMAP7 


Same as AGN, but with a WMAP7 cosmology 


DBLIMFV1618 


Top-heavy IMF at high pressure, extra SN energy in wind velocity 


DMONLY 


No baryons, cold dark matter only 


DMONL Y_ WMAP7 


Same as DMONLY, but with a WMAP7 cosmology 


MILL 


Millennium simulation cosmology (i.e. WMAPl), rj = 4 (twice the SN energy of REF) 


NOSN 


No SN energy feedback 


NOSN.NOZCOOL 


No SN energy feedback and cooling assumes primordial abundances 


NOZCOOL 


Cooling assumes primordial abundances 


WDENS 


Wind mass loading and velocity depend on gas density (SN energy as REF) 


WML1V848 


Wind mass loading rj = I, velocity v„ = 848 kms"-' (SN energy as REF) 


WML4 


Wind mass loading rj = 4 (twice the SN energy of REF) 



in Table [T] For more d etails and other models we again refer 
to lSchave et"aD l|2010h . 

The simulation DMONLY includes only dark matter, 
hence the only active physical process is gravity. This model 
is useful, as many (semi-) analytical models for the matter 
power spectrum assume that baryons are unimportant on 
large scales. 

The NOSN simulation excludes supernova feedback, 
and the simulation NOZCOOL assumes primordial abun- 
dances when computing cooling rates. The simulation 
NOSN.NOZCOOL excludes both SN feedback and metal- 
line cooling. Naturally, none of the three simulations can be 
considered realistic as we know that the omitted processes 
exist, but they are valuable tools to investigate on what 
scales and in what measure these processes affect the total 
matter power spectrum. In fact, the same may be said for 
the other models we consider as all, except for AGN, suf- 
fer from the overcooling problem and hence apparently miss 
an important process that does occur in nature (be it AGN 
feedback or something else). 

Supernova feedback models suffer from large uncertain- 
ties due to the limited resolution of the simulations and 
a lack of observational constraints. Though the product of 
the initial wind mass loading and the initial wind velocity 
squared, t]v^ , determines the energy injected into the winds 
per unit stellar mass and is therefore limited from above by 
the energy available from the SNe, the individual parameters 
are poorly constrained and can thus be varied. One varia- 
tion on the reference model that uses the same SN energy 
per unit stellar mass as REF is WMLl V848, in which the 
wind mass loading is reduced by a factor of 2 while the wind 
velocity is increased by a factor of \/2- Another such vari- 
ation is WDENS, in which the wind parameters scale with 
the density of the gas from which the star particle formed: 
the wind velocity as Ww oc nj/^, and the wind mass loading 
as oc v^'^ oc n^j^^'^. Both parameters are equal to their 
fiducial values for stars formed at the density threshold for 
star formation. For the polytropic EoS that we impose onto 
the ISM, the wind velocity in this model scales with the lo- 



cal effective sound speed, as might be the case for thermally 
driven winds. 

We also compare to models where the SN energy is var- 
ied. One scenario in which the SN energy may be higher 
than that in the reference model is when, under certain cir- 
cumstances, the IMF becomes top-heavy, meaning that rela- 
tively more high-mass stars are produced. It is expected that 
the IMF is top-he avy at high redshift and low metallicity 
(e.g. lLarsonlll998l ). and both observations and theory sug- 
gest that it may be top-heavy in extreme envi ronments like 
the galactic centre and s tarburst galaxies (e.g. iBaugh et al.l 
l200£l : lBartko et al.ll201ol ). In the simulation DBLIMFV1618, 
the latter effect is modelled by a switch from the Chabrier 
IMF to one that follows 4>{rn) oc once the gas reaches 
a certain pressure threshold, which is set so that ~ 10% of 
the stellar mass forms with a top-heavy IMF. In this case, 
the emissivity in ionizing photons goes up by a factor 7.3, 
and it is assumed that the SN energy scales up by the same 
factor. In the model we consider here, this extra energy is 
used to raise the wind velocity by a factor \/7.3. 

The final model that we consider that only differs from 
REF in terms of its wind parameters is WML4, in which 
the SN energy per unit stellar mass is doubled by simply 
increasing the wind mass loading by a factor of two. The 
same is done in the simulation MILL. However, the most 
important feature of the latter is that it uses the same val- 
ues for t he cosmological para meters as the Millennium sim- 
ulation (|Springel et al.ll2005l ). These are derived from first- 
year WMAP data and are given by: {ilm, fib, Oa, cts, '^s, h} 
= {0.25, 0.045, 0.75, 0.9, 1.0, 0.73}. 

The last and, for our purposes, most important physics 
variation we consider here adds a phenomenon that has 
proved to be increasingly necessary to reconcile theory and 
observations, from the scales of individual galaxies to clus- 
ters: Active Galactic Nuclei, or AGN. They are caused by 
the emission of large amounts of energy from the accret- 
ing supermassive black holes that reside at the centres of 
galaxies, in the form of radiation that may couple to the 
gas and relativistic jets caused by the magnetic field of the 
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infalling material, which can heat and displace gas out to 
very large distances. AGN have been invoked to explain, for 
example, the low star formation rates of high-mass galaxies 
and the suppression of cooling flows in clusters. Moreover, 
iLevine fc GnedinI ('2006') have used a toy model to demon- 
strate that AGN feedback may provide sufficient energy to 
have a large effect on the matter power spectrum. 

We model the growth of supermassive black holes and 
the assoc iated feedback processe s using the prescription de- 
tailed in 'Booth & Schaye ('2009"), which is an extension of 
that by[Springcl ct_al, ( 2005). During the simulation, a black 
hole seed particle with mass msecd = 9 X 10*/i-^Mq (i.e. 
10~^ m-baryon) is placed at the centre of every dark mat- 
ter halo whose mass exceeds mhaio.min = 4 x 10^" Mq 
(corresponding to 10^ dark matter particles). These parti- 
cles then accumulate mass from the surrounding gas at an 
(Eddingtor i-limited) rate based on Bondi-Hoyle-Lyttletq n 
accretion (|Bondi fc Hovld Il944l : iHovle fc LvttletonI Il939f) . 
but scaled up by a factor a to account for the lack of a 
cold gas phase and the finite numerical resolution. However, 
for densities below our star formation threshold we do not 
expect a cold phase to be present and we therefore set a 
equal to unity. To ensure a smooth transition, a is made to 
depend on the density of the gas: 

r 1 if riH < 

«=| (^i)" otherwise. 

Here the density threshold njj is the critical value required 
for the formation of a cold interstellar gas phase (rin = 
0.1 cm~^; see i]2.H) . Models of this type are called 'constant- 
P models', and the fiducial value /3 = 2 is used throughout 
this paper. 

The black holes inject 1.5 per cent of the rest mass en- 
ergy of the accreted gas into the surrounding matter in the 
form of heat. This feedback efficiency determines the nor- 
malisation, but not the slope, of th e relations betwee n black 
hol e mass and galaxy pro perties. iBooth fc Schavj (|2009l ) 
and lBooth fc Schavd (|201ll ') demonstrate that this efficiency 
reproduces the observed relations between BH mass and 
both stellar m ass and stellar ve locity dispersion, as well as 
their evolution. iMcCarthv et all ((2010) have shown that the 
AGN simulation, but not the reference model, provides ex- 
cellent agreement with both optical and X-ray observables 
of groups of galaxies at redshift zero. In particular, it repro- 
duces the temperature, entropy, and metallicity profiles of 
the gas, the stellar masses, star formation rates, and age dis- 
tributions of the central galaxies, and the relations between 
X-ray luminosity and both temperature and mass. We there- 
fore consider simulation AGN to be more realistic than our 
other models. As we shall see in ij3l the inclusion of AGN 
feedback greatly affects the power spectrum on a large range 
of scales. 

Finally, we have re-run two simulations, DM ONLY 
and AGN, with cosmo logical parameters d erived from the 
WMAP 7-year results (|Komatsu et al.ll201ll '): {Q^, fib, !^a, 
as, Wb, h} = {0.272, 0.0455, 0.728, 0.81, 0.967, 0.704}. These 
versions are called DM0NLY.WMAP7 and AGN.WMAP7, 
respectively. We consider the latter to be our most realis- 
tic and up-to-date model. Note that the linear input power 
spectra used for the initial conditions of these simulations 
have not been generated by CMBFAST, but by the more up- 



to-date f90 package camb (|Lewis fc Challinoill2002l . version 
January 2010). 

2.3 Power spectrum calculation 

The distribution of matter in the Universe can be described 
by a continuous density function, p(x), where the vector x 
specifies the position relative to some arbitrary origin. Given 
this density field, we consider fluctuations, 5{x), defined as: 

,(x).^Wzi, (2) 
P 

where p is the global mean density. We can relate this density 
contrast to wave modes 5k via a discrete Fourier transform: 

5(x) = ^4e-*''. (3) 

k 

The density field can thus be seen as made up of waves with 
certain amplitudes and phases, with wave vectors k. We now 
define the power spectrum, P{k), as: 

P{k)^v{\5^\')^, (4) 

where V is the volume under consideration. The power spec- 
trum is therefore obtained by collecting the amplitudes- 
squared of all wave modes with the same wave number 
k = |k|, and averaging them. This makes it clear that the 
power spectrum is a statistical tool, whose accuracy in- 
creases when more waves of the same length are available 
(i.e. when the scale 2n/k is small compared to the size of 
the box). We will present our results using what is often 
called the dimensionless matter power spectrum, which is 
defined as: 

A'(fc) = ^Pik). (5) 

The dimensionless power spectrum scales with the mass vari- 
ance, cr^(M), where M = ^ (^)^ P- Note that using A^(/c) 
instead of P{k) does not affect the relative differences be- 
tween power spectra. 

The code we have chosen to use to obtain accurate 
power spectrum estimations from our simulat ions is the pub- 
licly available f90 package called powmes (jColombi et al.l 
boogf ). The advantages of powmes stem from the use of the 
Fourier- Taylor transform, which allows analytical control of 
the biases introduced, and the use of foldings of the particle 
distribution, which allow the dynamic range to be extended 
to arbitrarily high wave numbers while keeping the statisti- 
cal errors b ounded. For a full des cription of these methods 
we refer to IColombi et all l|2009f ). We have compared the 
performance of powmes with respect to power spectrum 
estimators using simple NGP, CIC and TSC interpolation 
schemes, and found that powmes is capable of obtaining 
far more accurate power spectra over a larger spectral range 
within the same computation time. 

We have expanded powmes with the possibility to con- 
sider only one group of particles at a time, in order to see 
which parts of the power spectrum are dominated by the 
contribution of, for example, cold dark matter (see ijS.Sp . 
Finally, we performed extensive timing tests using different 
grids, foldings, CPUs and particle numbers which, combined 
with the performance results from lColombi et al] (|2009l ). re- 
sulted in the fiducial values G — 256 and F = 7 for the 
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number of grid points on a side and the number of foldings, 
respectively. 

2.3.1 Discreteness and other numerical limitations 

In Appendix A we demonstrate that the simulations are suf- 
ficiently converged with respect to increases in the numerical 
resolution to predict the power spectrum with better than 
1% accuracy for k < lO/iMpc"^. This range is greatly ex- 
panded in both directions if we only consider the relative 
differences in power between simulations. 

Besides numerical resolution, the predicted power spec- 
tra may be affected by sample variance, which is generally 
called cosmic variance in cosmology. This is caused by the 
finite volume of the box and by the fact that each simulation 
provides only a single realisation of the underlying statisti- 
cal distribution. Note that finite volume effects are different 
for the simulations than for observational surveys, because 
the mean density in the simulation boxes is always equal to 
the cosmic mean. In Appendix A we show that finite volume 
effects may cause us to underestimate the effects of baryons 
on scales of several tens of Mpc, i.e. close to the size of the 
box. While the fact that we only use a single realisation of 
the initial conditions prevents us from obtaining highly ac- 
curate absolute values for the power spectrum on scales close 
to the box size, it does not prevent us from investigating the 
relative changes in power caused by baryon physics. 

Finally, we are limited in our determination of the power 
spectrum by the discreteness of the density field. Because we 
use particles to represent a continuous field, there will always 
be non-zero power present at all scales, called white noise 
or shot noise. If we assume the particle distribution to be 
a local Poisson realisation of a stationary random field, an 
assumption used in any calculation of the power spectrum 
and one that is expected to be valid for an evolved distri- 
bution0 this white noise component can be calculate d (see, 
for example, Peebles 1980, 1993; Co lombi et al.ll2009l '). Sub- 
tracting the shot noise from the initial estimate of the power 
spectrum will make the latter somewhat more accurate, but 
one should still expect the uncertainty on the estimate of the 
power spectrum to increase dramatically when the intrinsic 
power spectrum falls far below the shot noise level. The con- 
tribution of shot noise to P{k) is independent of k. Hence, 
if we use A^(fe) as the measure of the power spectrum, then 
the shot noise level will scale as k^. In the following sec- 
tion, the scale at which the shot noise of each simulation is 
equal to (the white noise corrected) A^(fe) is denoted by a 
circle, while it is shown explicitly in Appendix A. Note that 
the theoretical shot noise level has been subtracted for all 
power spectra shown in this paper. 
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Figure 1. Comparison of the matter po wer spectrum of 
DMONLY.L100 N512 with an alytical fits by IPeacock &: Doddd 
(1996, PD96) and lSmith et al.l ll2003. HALOFIT) at redshift zero. 
The small circle, drawn in this and all following plots showing 
A^(fc), indicates the scale below which the (subtracted) shot noise 
in the simulation becomes significant, and the dashed purple curve 
shows the linear input power spectrum of the simulations. The 
bottom panel shows the ratios of the power spectra from theoret- 
ical models and the simulation. There is good agreement down 
to scales of a few Mpc, especially for the more recent HALOFIT 
model, but on smaller scales DMONLY predicts up to twice as 
much power as HALOFIT. For A < 10^ h'^ kpc the power in the 
DMONLY simulation drops due to a lack of resolution. 



of our dark matter only simulation to predictions from the 
literature. We investigate the effects of adding or modify- 
ing prescriptions for baryonic processes in i]3.2l We examine 
how well CDM, gas, and stars trace each other and con- 
sider the contributions of these different components to the 
total power in the reference simulation in i]3.3l and we ex- 
amine the back-reaction of baryons on the CDM for the two 
most important simulations, REF and AGN, in ^3.4\ Fi- 
nally, in ^3.5\ we take a closer look at model AGN, which 
we consider to be our most realistic simulation because it 
reproduc es the optical and X-ra y observations of groups of 
galaxies (jMcCarthv et al.1 l2010l ) . We investigate the effect 
of using the WMAP7 rather than the WMAP3 cosmology, 
compare to widely used model power spectra, and consider 
the evolution of the effect of baryons on the matter power 
spectrum. 



3 RESULTS 

In this section we present the power spectra obtained from 
our simulations. In i]3.1l we compare the power spectrum 



^ The discreteness noise can initially be much smaller if the par- 
ticles are arranged on a grid or in a "glass-like" fashion. Particles 
in low-density regions may retain memory of their initial distri- 
bution, reducing the noise below the level expected for a Poisson 
distribution. 



3.1 Comparison of a dark matter only simulation 
to models 

In this section we compare the power spectrum of our 
DMONLY simulation to thos e predi cted by the widely used 
models of iPeacock fc Doddd l|l996l . hereafter PD96) and 
ISmith et al l l|2003l . hereafter HALOFIT) . 

The PD96 model is an extension of what is known as the 
HKLM model (Hamilton ct al. 1991), which first introduced 
a universal analytical formula to map the linear correlation 
function into a non-linear one, the coefficients of which were 
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estimated using N-body simulations. Both of these models 
assume spherical collapse of fluctuations that have reached 
a certain overdensity, followed by stable clustering (which 
states that the mean physical sep aration of particles is con - 
stant on sufficiently small scalesl . [Peacock fc Doddi l|l994h . 
followed by PD96, expanded on the groundwork laid by 
HKLM by presenting a version of the method that worked 
with power spectra instead and allowed for Q ^ l,a non-zero 
cosmological constant and large negative spectral indices. 

However, numerical simulations have shown that the as- 
sumption of stable clusteri ng is not always va lid. The more 
recent HALOFIT model by|SmitheiaD ((2OO3I ') aimed to im- 
prove on PD96 by taking this into account. This method is 
based on concepts from the "halo model" , in which the den- 
sity field is viewed as a distribution of isolated haloes (e.g . 
IPeacock fc Smithll2000l : ISeliakll2000l : ICoorav fc Shethll2002h . 
It is then assumed that the power spectrum can be split into 
two parts: a large-scale quasi-linear term that is due to the 
clustering of separate haloes (the 2-halo term), and a small- 
scale term caused by the correlation of subhaloes within the 
same parent halo (the 1-halo term). Their resulting analyti- 
cal formulae were fit to power spectra obtained from N-body 
simulations. 

To create power spectra that conform with these mod- 
els and the cosmological parameters used in our simulations, 
we have u tilised the publi c ly ava ilable package iCoSMO, de- 
scribed in lRefregier et al.1 (|2008|'). We chose to generat e the 
linear power spectra using the Isisenstein fc Hul (Il999l. EH) 
transf er function. We have also tried using the' Bardeen et al] 
l|l986l . BBKS) transfer function to generate initial conditions 
for the PD96 model, as this is the one originally used by the 
authors, which introduced only minor differences with re- 
spect to the results shown here (1 — 10% lower power for 
k > lO/iMpc"^). 

In Figure [T] we compare these models to the simula- 
tion that, like the theoretical models for non-linear growth, 
includes only dark matter (DMONLY). For reference, the 
dashed curve shows the linear input power spectrum of the 
simulations. The bottom panel shows the ratio of the ana- 
lytical predictions to our results. Note that we have omit- 
ted the first wave mode (at A = 100 Mpc) in all of our 
figures because we cannot sample the power spectrum on 
the scale of the simulation box. We see that the dark mat- 
ter power spectrum follows the analytical predictions pretty 
well on large scales (except on the scale of the simulation 
box), and that HALOFIT provides a better match than the 
PD96 model, as expected. However, on scales below a few 
Mpc the theoretical models start to severely underestimate 
the amount of structure formed in the simulation, and the 
difference between HALOFIT and the DMONLY simula- 
tion reaches a factor of 2 on scales of 1 — 3 x 10~^ h'^ Mpc. 
The rapid decline of the DMONLY power spectrum for 
k > 100/iMpc"^ (A < 10^^"^kpc) is mostly due to the 
underproduction of low-mass haloes due to the finite resolu- 
tion (see Appendix A). While we will always show the power 
spectrum up to « 500/iMpc~^, we are mainly interested 
in the scales relevant for upcoming surveys, k < lOh Mpc~^. 
As discussed in Appendix A, for k ^ 10/iMpc~^ numerical 
convergence may become an issue. Note that the power spec- 
trum of the simulation remains reasonably well-behaved far 
below the theoretical shot noise level (i.e. well to the right 
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Figure 2. A comparison of the total matter power spectra 
of DMONLY.L100N512 (black), REF.L100N512 (green) and 
AGN.L100N512 (red), at redshift 2 = 0. The bottom panel 
shows the absolute value of the relative difference of the latter 
two with respect to DMONLY; solid (dashed) curves indicate 
that the power is higher (lower) than for DMONLY. The dotted, 
horizontal line shows the 1% level. Note that the first wave mode 
has been omitted as it holds no information. While pressure forces 
smooth the baryonic density field on intermediate scales, cooling 
allows the baryons to increase the total power on small scales. 
The addition of AGN feedback, which is required to match ob- 
servations of groups, has an enormous effect, reducing the power 
by > 10% for k > 1/iMpc-i. 



of the small circle), indicating that the subtraction of this 
noise component is fairly accurate. 

Newer implementations of the halo model exist, based 
on fits to more recent N-body simulations. These mod- 
els improve on the HALOFIT model by including a vari- 
able concentration-mass relation ( such as those derived by 
iNeto et al.ll2007l : iDuffv et al1l2008h and have been shown to 
reproduce the power spectr a from simulations with higher 
accuracy (e.g. iHilbert et all [2009,). Since no suitable codes 
employing these models were ava ilable, we do not com pare 
to their results here. However, as iHilbert et af] (|2009l ) have 
shown that usin g the halo model with the concentration- 
mass relation of iNeto et al.l (|2007l ) increases the power at 
intermediate scales, we suspect that such models would pro- 
vide a better match to the power spectrum of DMONLY. 

3.2 The relative effects of different baryonic 
processes 

In this section we present our main results, demonstrating 
how single baryonic processes, or implementations thereof, 
can influence the matter power spectrum. While we will fo- 
cus mainly on the range of scales relevant to upcoming weak 
lensing surveys, 0.1/iMpc~^ <k < lO/iMpc"^, we will also 
discuss the differences at the much smaller scales that our 
simulations allow us to probe. We note again that all power 
spectra are taken from simulations with L = 100 Mpc 
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and TV = 512 at redshift zero, and that, unless stated other- 
wise, all simulations are evolved from the same initial con- 
ditions. 

We start by comparing the power spectra of DMONLY, 
the reference simulation (REF) and AGN, in Figure [2] The 
panel at the bottom of most plots in this section shows 
the absolute value of the relative difference between power 
spectra. The dotted, horizontal line shows the 1% level: any 
differences above this level will thus affect the statistics of 
surveys that aim to measure the power spectrum to this 
accuracy. 

It is immediately clear from the comparison between 
DMONLY and REF that the contribution of the baryons 
is significant, decreasing the power by more than 1% for 
k « 0.8 — 5/iMpc~^. This is because gas pressure smooths 
the density field relative to that expected from dark matter 
alone. On scales smaller than 1 /i"^ Mpc [k > 6/iMpc~^), 
the power in the reference simulation quickly rises far above 
that of the dark matter only simulation, because radiative 
cooling enables gas to cluster on smaller scales than the dark 
matter. These results confirm the findings of previous stud- 
ies, a t least qualitatively fe.g. [jing et al]l2006l : iRudd et al.l 
l2008l : iGuillet et aLll201Cll : ICasarini et al."201l'). 

However, when AGN feedback is included, the results 
change drastically. In this case, the reduction in power rel- 
ative to DMONLY already reaches 1% for k ^ 0.3 /iMpc"^ 
(A '^2Qh~^ Mpc) and exceeds 10% for 2 < fc < 50 /i Mpc"^ 
We thus see that AGN feedback even suppresses the to- 
tal matter power spectrum on very large scales. The enor- 
mous effect of AGN feedback is due to the removal of gas 
from (groups of) galaxies. That large amounts of gas are in- 
deed being moved to la rge radii in this sim ulation has been 
shown by, for exam ple, Duf fv et al.l (|201G| . e.g. Figures 1 & 
2) and iMcCarthv et al. ( 201ll . e.g. Figure 3). Because the 
AGN reside in massive and thus strongly clustered objects, 
the power is suppressed out to scales that exceed the scale 
on which individual objects move the gas. 

Figure |3] shows the difference in the power spectra pre- 
dicted by a variety of simulations relative to that predicted 
by the reference simulation. The models are listed in Table[T] 
and were described in ii2.2l The top panel shows the effect 
of turning off SN feedback and/or metal- line cooling. Since 
SN feedback heats and ejects gas, we expect it to decrease 
the small-scale power. Indeed, the power in NOSN is > 1% 
higher than in the reference simulation for k > 4/iMpc~^ 
and the difference reaches 10% at fc ~ lO/iMpc"^. The ab- 
sence of SN feedback also increases the star formation rate, 
making stars the dominant contributor to the total matter 
power spectrum out to larger scales (not shown). 

Turning off metal-line cooling reduces the power on 
small scales because less gas is able to cool down and accrete 
onto galaxies. Indeed, model NOZCOOL predicts 10 — 50% 
less power for k > 30/iMpc~^. However, the absence of 
metal-line cooling increases the power by several percent for 
A ~ l/i-^Mpc because the lower cooling rates force more 
gas to remain at large distances from the halo centres. 

Even though the effects of SN feedback and metal-line 
cooling are somewhat opposite in nature, as the former in- 
creases the energy of the gas while the latter allows the gas to 
radiate more of its thermal energy away, removing both pro- 
cesses in the simulation NOSN^NOZCOOL still introduces 
differences of about 1 — 10% for k > 2/iMpc~^ relative to 
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Figure 3. Comparisons of z = power spectra predicted by simu- 
lations incorporating different physical processes to that predicted 
by the reference simulation. The panels are similar to the bottom 
panel of Figure [2] but now show differences relative to REF. The 
thin black curve that is repeated in all panels shows the relative 
difference with DMONLY. Colours indicate different simulations, 
while different line styles indicate whether the power is reduced 
or increased relative to the reference simulation. 
Top: A simulation without SN feedback (blue), one without 
metal-line cooling (green) and one that excludes both effects 
(red). SN feedback decreases the power on all scales. Metal-line 
cooling decreases the power for A > 0.4 Mpc but increases the 
power on smaller scales. The effects of removing both SN feed- 
back and metal-line cooling are > 10% for k > 20/iMpc~^ and 
> 1% for k > 2hMpc-i. 

Middle: Different SN wind models which all use the same amount 
of SN energy per unit stellar mass (see text). The effects of vary- 
ing the implementation of SN feedback, while keeping the SN 
energy that is injected per unit stellar mass the same, are > 10% 
for k > 10/iMpc-l and > 1% for k > IkMpc'^. 
Bottom: Models with different feedback energies and processes, 
see text for details. Including a top-heavy IMF at high pres- 
sure {DBLIMFV1618) or AGN feedback (AGN) greatly reduces 
the power. The reduction caused by the latter is > 10% for 
k > 2hMpc-i and > 1% for k > 0.4/iMpc-l. 
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Figure 4. Difference of the 2 = matter power spectrum in a 
simulation using a WMAPl cosmology (MILL) relative to that 
of the REF model, which assumes the WMAP3 cosmology, af- 
ter rescaling the former to match the latter on the scale of the 
simulation box (A = 100/i~^Mpc, not shown). WML4 is shown 
for reference as this simulation uses the same baryonic physics as 
MILL. For k>3h Mpc~^, the effect of AGN feedback is at least 
as strong as that of this unrealistically large change in cosmology. 



the reference simulation. It is therefore vital to take both SN 
feedback and metal-line cooling into account if one wants to 
predict the matter power spectrum with an accuracy better 
than 10%. 

We compare models that use different prescriptions for 
SN feedback, but the same amount of SN energy per unit 
stellar mass as REF, in the middle panel of Figure |31 In 
WMLl V848 the SN energy is distributed over half as much 
gas, but the initial wind velocity is a factor \/2 higher, re- 
sulting in more effective SN feedback in all but the low- 
est mass galaxies. The differences with respect to the ref- 
erence model extend to even larger scales than when SN 
feedback is removed entirely: the power is reduced by > 1% 
for k > 1/iMpc-i and by > 10% for k > lOfeMpc^^ In 
model WDENS the initial wind velocity increases with the 
local sound speed in the ISM, but the mass loading is ad- 
justed so as to keep the amount of SN energy per unit stellar 
mass equal to that in REF. This implementation results in 
an even stronger decrease in power on scales < 10^"^ Mpc. 
In both these models, the reduction in power is caused by 
the increased effectiveness of SN feedback in driving outflows 
of gas. We stress that because of our lack of understanding 
of the effects of SN feedback, there is a priori no reason to 
assume that the model used in the reference simulation is a 
better approximation to reality than the models we compare 
to here. 

In fact, it is possible that the SN energy per unit stellar 
mass is different from the value assumed in the REF model, 
or that it varies with environment. Model DBLIMFV1618, 
which we compare with REF in the bottom panel of Fig- 
ure [3l uses a top-heavy IMF in high-pressure environments. 
Such an IMF yields more SNe per unit stellar mass which 
decreases the power by > 1% for k > 0.7/iMpc"^ and by 
> 10% for k > 4 ft Mpc Clearly, it will be necessary to 
understand any environmental dependence of the IMF in or- 



der to predict the matter power spectrum to 1% accuracy 
on the scales relevant for upcoming surveys. 

On the other hand, doubling the wind mass loading, 
while keeping the wind velocity fixed to the value used in 
REF, as is done in WML4, has a far more modest effect. 
This is because the wind velocities are too low to signif- 
icantly disturb the high-pressure ISM of massive galaxies. 
The differences with respect to the reference model are lim- 
ited to < 1% for k < 10/iMpc-\ 

The bottom panel of Figure [3] also compares the ref- 
erence simulation to model AGN, which differs from REF 
by the inclusion of a phenomenon that has been shown 
to play a role in many contexts and that strongly im- 
proves the agreement w ith observations of groups of galaxies 
(iMcCarthv et al.ll2010l ). Like SN feedback, AGN feedback 
decreases the power by heating and ejecting gas, but the 
effect is more dramatic than that of the standard SN feed- 
back model, both in scope and magnitude. With respect to 
the reference model, the power is decreased by > 30% for 
k > 10/iMpc"^ and by > 5% for k > l/iMpc'\ The re- 
duction in power only falls below 1% for k < 0.4/iMpc~^ 
(A > 10 h'^ Mpc). Note that the effect of AGN feedback is 
strikingly similar to, albeit stronger than, that of the stellar 
feedback model that uses a top-heavy IMF in high-pressure 
environments. 

It is clear that many different baryonic processes, and 
even slightly different implementations thereof, are capable 
of introducing significant differences in the matter power 
spectrum on scales relevant for observational cosmology. To 
put the effects of baryons into perspective, we compare to 
a simulation with a very different cosmology, MILL, in Fig- 
ure [31 The difference between the cosmology derived from 
the first-year WMAP data used in MILL and the one derived 
from the 3-year WMAP data used in the other simulations 
is large; in fact, the difference is much larger than the er- 
ror bars of the most recent data allow. For refere nce, we 
note that the currently favoured cosmology (Komats u et al.l 
l201lf ) lies in between those given by WMAPl and WMAP3. 
To account for the difference in normalisation of the MILL 
power spectrum, which is caused mainly by its higher Qm 
and (Tg values, we have rescaled it to have the same power at 
the box size as REF. Still, the effect on the power spectrum 
exceeds 10% for k > 0.2ftMpc"^ A quick comparison with 
WML4, which uses the exact same baryon physics as MILL 
and twice the SN wind mass loading used in REF, shows 
that the effect of the change in mass loading is relatively 
small, as we had already shown in Figure [31 However, we 
see that for k >3h Mpc~^, the effect of AGN feedback is at 
least as strong as that of this unrealistically large change in 
cosmology. We thus conclude that baryonic effects are not 
only significant at the ~ 1% level, but can even be larger 
than a "very wrong" choice of cosmology. 

Almost all theoretical models used in the literature con- 
sider only CDM, assuming that the baryons follow the dark 
matter perfectly for k < l/iMpc~^. We have shown (see 
Fig. [21) that the fact that baryons experience gas pressure 
reduces the power on large scales, while their ability to radi- 
ate away their thermal energy increases the power on small 
scales. If we ignore AGN feedback, as has been done in 
all previous work, we find that the power is reduced by at 
least a few percent for 0.8 < k < 5/iMpc~^ and that the 
power is increased for k > 7 ft Mpc" ^, with the difference 
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Figure 6. The back-reaction of baryons on the CDM. The blue 
curves show the relative difference between the power spectrum 
of the CDM component, after scaling the CDM density by the 
factor $7m/(f^m — f^b)i and that of a dark matter only simu- 
lation for either the REF (top panel) or AGN (bottom panel) 
model. For comparison, the relative differences between the total 
matter power spectra of the baryonic simulations and DMONLY 
is shown by the black curves. Baryons increase the small-scale 
power in the CDM component. However, when AGN feedback is 
included, the power in the CDM component drops 1 — 10% be- 
low that of the DMONLY simulations for 0.2h-lMpc < A < 
2h~^ Mpc. 



reaching approximately 6% at fc = 10/iMpc^^ for the refer- 
ence modeL However, the single process of AGN feedback, 
which improves the agreement with observations of groups of 
galaxies, reduces the power by > 10% over the whole range 
1 < fc < 10 /iMpc~^ and the reduction only drops below 1% 
for k < 0.3 /iMpc~^. Highly efficient SN feedback, as may for 
example result from a top-heavy IMF in starbursts, would 
have nearly as large an effect. One can therefore not expect 
to constrain the primordial power spectrum more accurately 
until such processes are better understood and included in 
theoretical models. 



3.3 Contributions of dark matter, gas and stars 

Generally, power spectra are calculated using all matter in- 
side the computational volume. This total matter power 
spectrum is what is measurable using e.g. gravitational lens- 
ing surveys. However, as we have a larger freedom of mea- 
surement using simulations, we can also consider the power 
in different components, for example to see which parts of 



the power spectrum are dominated by baryonic matter or 
how baryons change the distribution of cold dark matter. 

On sufficiently large scales the baryons will trace the 
dark matter. Hence, when averaged over these scales, the 
baryonic and CDM densities are given by 

Pcdm ~ Ptot, 

Pbar = TT—Ptot- (O) 

We can now use these expressions to estimate the relative 
contributions of correlations between particle types to the 
total matter power spectrum. Using Ptot(fc) (|ptot(fc)P) oc 

(|/5cdm|^) + (PcdmPbar) + (PcdmPbar) + (|/5bar|^), WC find, for 

sufficiently small k: 

Pec = TT^ ^Ptot ~ 0.68Ptot , 

Pcb + Pbc = ^^^^^%^^Ptot^0.29Ptot, (7) 

Pbb = T^Ptot ~0.03Ptot. 

Hence, on large scales we expect the power due to the auto- 
correlation of CDM to dominate the total matter power 
spectrum, with a significant contribution from the cross 
terms Pcb and Pbc- 

The four panels of Figure [5] show power spectra for the 
REF.L100N512 (left) and AGN.L100N512 (right) simula- 
tion at 2 = 0, both for the total matter (solid black) and 
for individual components (coloured curves). For reference, 
we also show the power spectrum for DMONLY_L100N512 
(dashed black). The top row shows the power spectra of 
5i = {pi — pi) I pi. This definition ensures that the power 
spectra of all components i converge on large scales, which 
allows us to examine how well different components trace 
each other. The bottom row, on the other hand, shows the 
power spectra of (5^ = {pi — ptot)//Otot, which allows us to 
estimate the contributions of different components to the 
total matter power spectrum. 

Looking at the top-left panel, we see that, as expected, 
the baryonic components trace the dark matter well at 
the largest scales. However, significant differences exist for 
A < 10 hr^ Mpc. Observe that, at scales of several hundred 
kpc and smaller, the difference between REF and the dark 
matter only simulation is larger than that between the lat- 
ter and the analytical models we compared to earlier (see 
Fig. [T]). In fact, the difference between the cold dark matter 
component of the reference simulation and DMONLY is also 
larger than that between the latter and the analytic models. 
This is due to the back-reaction of the baryons on the dark 
matter, which we will discuss in i ]3.4l 

Next, we turn to the bottom-left panel of Figure[5]which 
shows that cold dark matter dominates the power spectrum 
on large scales, as expected, although the contribution from 
the CDM-baryon cross power spectrum (not shown) is im- 
portant as well. The contribution of baryons is significant 
for A < 10^ kpc and dominates below 60/i~^kpc. The 
strong small-scale baryonic clustering is the direct conse- 
quence of gas cooling and galaxy formation. Taking a look 
at how the baryonic component is itself built up, we see 
that gas dominates the baryonic power spectrum on large 
scales, but that stars take over for A < lh~^ Mpc. The gas 
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Figure 5. Decomposing the z = total power spectra (black) into the contributions from cold dark matter (blue), gas (green) and 
stars/black holes (red). The left and right columns show results for REF^L100N512 and AGN^L100N512. In the top row the density 
contrast of each component i is defined relative to its own mean density, i.e. 5i = (pi — pi)/pi. This guarantees that all power spectra 
converge on large scales, thus enabling a straightforward comparison of their shapes. In the bottom row the density contrast of each 
component is defined relative to the total mean density, i.e. Si = (pi — ptot)/ptot, which allows one to compare their contributions to the 
total power. The power spectrum of the gas flattens or even decreases for A < 1 h"^ Mpc as a result of pressure smoothing, but its ability 
to cool allows it to increase again on galaxy scales (A < 10^ kpc). The power spectrum of the stellar component, which is a product 
of the collapse of cooling gas, increases most rapidly towards smaller scales. While stars dominate the total power for A <C 10^ kpc 
in REF, dark matter dominates on all scales when AGN feedback is included. 



power spectrum flattens for A < l/i~^Mpc, which corre- 
sponds to the virial radii of groups of galaxies, but steepens 
again for A < 0.1 Mpc, i.e. galaxy scales. The reason for 
the decrease in slope around lh~^ Mpc is threefold. First, 
the pressure of the hot gas smooths its distribution on the 
scales of groups and clusters of galaxies. Second, as the gas 
collapses it fragments and forms stars. Third, due to stel- 
lar feedback the gas is distributed out to large distances, 
reducing the power. 

The inclusion of AGN feedback greatly impacts the 
matter power spectrum on a wide range of scales. Compar- 
ing the top panels of Figure [S] we see that AGN feedback 
strongly decreases the power in the gas and stellar compo- 
nents relative to that of the dark matter for A < lh~^ Mpc. 
A comparison of the bottom panels reveals that the contri- 
bution of stars to the total power is reduced the most, with 
the reduction factor increasing from an order of magnitude 
on the largest scales to more than two orders of magnitude 
on the smallest scales. This clearly shows that AGN feed- 
back suppresses star formation, as required to solve the over- 
cooling problem. For the gas component the change is also 
dramatic. While I^l^^{k) = 1 for A ~ 3/i"^Mpc in REF, 



this level of gas power is only reached at 100 /i"^ kpc for 
AGN. The suppression of baryonic structure by AGN feed- 
back makes dark matter the dominant component of the 
power spectrum on all scales shown, although it is important 
to note that the dark matter distribution is also signiflcantly 
affected by the AGN, as we shall see next. 

3.4 The back-reaction of baryons on the dark 
matter 

Even though dark matter is unable to cool through the emis- 
sion of radiation, its distribution can still be altered by the 
inclusion of baryons due to changes in the gravitational po- 
tential. We examine this back-reaction of the baryons on 
the dark matter for the reference and AGN simulations in 
the left and right panels of Figure |S] respectively. In order 
to make a direct comparison, we have rescaled the density 
of the dark matter component of the simulations that in- 
clude baryons by multiplying it by the factor r2ni/(f2m — f^b)- 
The blue curve shows the relative differences between the 
power spectrum of the rescaled CDM component and that 
of DMONLY. 
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Figure 7. The dependence of the effect of AGN on cosmology. 
The curves show the relative differences between the z = matter 
power spectra for models AGN and DMONLY for our fiducial 
WMAP3 cosmology (green) and for the WMAP7 cosmology (red). 
Changing the cosmology has little impact on the relative effect of 
the baryonic processes. 
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Figure 8. Evolution of the relative difference between the matter 
power spectra of DM0NLY.WMAP7 and AGN.WMAP7. From 
red to blue, redshift decreases from 3 to zero. The erratic be- 
haviour of the 2 = 2 and 2 = 3 power spectra at the very smallest 
scales shown is due to a lack of resolution. For A > 1 Mpc the 
reduction in power due to baryons evolves only weakly for 2 < 1, 
but the transition from a decrease to an increase in power keeps 
moving to smaller scales. 



On scales k > 2/iMpc~^, corresponding to spatial 
scales A < 3/i~^ Mpc, the power in CDM structures in the 
reference simulation is increased by > 1% with respect to 
DMONLY. The difference continues to rise towards higher k, 
reaching 10% around k = 10 /i Mpc"'^. Because the baryons 
can cool, they are able to collapse to very high densities, and 
in the process they steepen the potential wells of virialized 
dark matter haloes, causing these to contract. The effect is 
larger closer to the centres of these haloes, i.e. on smaller 
scales. 

The back-reaction is quite different when AGN feed- 
back is included0 The dark matter haloes still contract on 
small scales, albeit by a smaller amount, but the power in 
the dark matter component of the AGN simulation is de- 
creased for scales > 200 kpc, corresponding to the sizes 
of haloes of L» galaxies. The reduction in the power of the 
CDM component in model AGN relative to DMONLY in- 
creases from roughly 1% at fc = 3/iMpc~^ to almost 10% 
around k = 10/iMpc~^. AGN-driven outflows redistribute 
gas to larger scales, which reduces the baryon fractions in 
haloes and results in shal lower potent i al we lls. This is con- 
sistent with the results of iDuffv et al.l (|2010l '). who used the 
same simulation to show that AGN feedback decreases the 
concentrations of dark matter haloes of groups and clusters. 
Note, however, that because AGN can drive gas beyond the 
virial radii of their host haloes, their effect on the power 
spectrum cannot be fully captured by a simple rescaling of 
the halo concentrations. 



° The small difference in power between the CDM component 
of AGN and DMONLY near the size of the box is most likely 
caused by errors in the power spectrum estimation. 



3.5 A closer look at the effects of AGN feedback 

In this section we examine our most realistic model for the 
baryonic physics, AGN, more closely. 

3.5.1 Dependence on cosmology 

Figure [7] shows how the relative difference between the 
z = power spectra of models AGN_WMAP7 and 
DMONLY. WMAP7, both of which use the WMAP7 cosmol- 
ogy, compares to that between the same physical models in 
the WMAP3 cosmology (the latter case was already shown 
in Figure [2J. Even though the power spectra are themselves 
strongly influenced by, for example, the much higher value of 
as in the WMAP7 cosmology, the relative change in power 
due to baryons is nearly identical, at least so long as AGN 
feedback is included. This is good news for observational 
cosmology. It means that, once the large current scatter in 
implementations of subgrid physics has converged, it may be 
possible to separate the baryonic effects from the cosmologi- 
cal ones when modelling the matter power spectrum. It also 
means that we can assume that our results of the previous 
sections, which were based on the WMAP3 version of the 
AGN simulation, apply also to model AGN.WMAP7. 

3.5.2 Evolution 

Next, we use the AGN.WMAP7 simulation, which we con- 
sider to be our most realistic model, to investigate the de- 
pendence of the effect of baryon physics on redshift. Figure[8] 
shows the relative difference between the power spectra of 
DM0NLY_WMAP7 and AGN.WMAP7 at redshifts 3, 2, 1, 
0.5 and zero. We see from this plot that on large scales, 
X > lh~^ Mpc, the reduction in power due to the gas does 
not evolve much for 2 < 1, although the differences between 



© 2011 RAS, MNRAS OOO.ITlfTSl 



Galaxy formation and the matter power spectrum 13 



the different redshifts remain large compared with the preci- 
sion of upcomi ng surveys. The w e ak ev olution below 2: = 2 is 
consistent with [McCarthy et alj (|201ll ). who found that the 
expulsion of gas due to AGN feedback takes place primar- 
ily at 2 < 2; < 4. On scales below 1 Mpc, on the other 
hand, the effects of baryonic processes on the power spec- 
trum keep increasing with time, with the transition point 
between a decrease and an increase in power steadily moving 
towards smaller scales. This is probably because the ejection 
of low-entropy halo gas at high redshift {z > 2) results in an 
increase of the entropy, and thus a r eduction of the coolin g 
rates, of hot halo gas at low redshift l|McCarthv et al.ll201lh . 



4 COMPARISON WITH PREVIOUS WORK 

Our predictions for the effect of baryons on the matter power 
spectrum agree qualitatively with those of other authors, 
provided we restrict ourselves to including the same baryonic 
feedback processes as were considered in those studies, flow- 
ever, previous simulations did not include AGN feedback 
and hence suffered from overcoolingl^l As we have demon- 
strated, AGN feedback (or very efficient stellar feedback) 
has a dramatic effect on the matter power spectrum over a 
large range of scales. In this section we will consider both the 
qualitative and quantitative differences with respect to pre- 
vio us wo rk, and examine how these may have come abo ut. 

[jiiig et all l\206(h used gadget ii (|Springell lioosi ') to 
run a simulation with a 100 ft"'^ Mpc box and 512^ gas 
and DM particles. Their simulation i ncluded radiative cool- 
ing an d star formation, and used the lSpringel fc Hernguistl 
l|2003h sub-grid model for the multiphase ISM and for galac- 
tic winds driven by star formation. Metal-line cooling and 
AGN feedback were not considered. They found that the 
power at fc = l/iMpc"'^ is reduced by ~ 1% relative to 
a dark matter only simulation at z = 0, which matches 
our results for the reference simulation very well. Further- 
more, in agreement with our reference model, they find that 
the inclusion of baryons increases the power by ~ 10% at 
k = 10/iMpc~^. However, they find that the transition from 
a relative decrease to a relative increase in power occurs at 
k~ 2h Mpc~^, while we find that it lies at k ~ 6h Mpc"'^. 

As the simulation of iJing et al.l (|2006D excludes metal- 
line cooling, we expect their results to be in better agreement 
with our own results for NOZCOOL. The main difference 
with respect to the reference simulation turns out to be the 
position of the transition point from a relative decrease to a 
relative increase in power, which shifts to fc « 2 — 3 /i Mpc~^ 
when metal-line cooling is turned off. Hence, using the sim- 
ulation NOZCOOL, w e reproduc e both the qualitative and 
quantitative results of ljing et al.l (|2006| ). even though bary- 
onic processes such as SN feedback are not implemented in 
the same way. 

iRudd et all l|2008h used the ART code l|Kravtsovlll999l ') 
ex panded with the Euler ian hydrodynamics solver described 
in lKravtsov et al.l l|2002l ). They used a 60 /i"^ Mpc box with 



The toy model of iLevine fc GnedinI l l2006h . which we briefly 
describe later in this section, did demonstrate, based purely on 
energetic grounds, that AGN feedback has the potential to have 
a large eflect on the matter power spectrum. 



256"^ particles, and included radiative cooling and heat- 
ing, metal-line cooling, star formation, therm al SN feedback 
(which is described in iKravtsov et al.l |2005| ) and chemical 
enrichment. AGN feedback was not considered. The effect 
of the baryons on the matter power spe ctrum they found 
is far more dramatic than that found bv ljing et al.l (|2006l ) 
and ourselves: a decrease in power of up to ~ 10% relative 
to a dark matter only simulation for k < l/iMpc~^, and a 
relative increase in power at fc > 1/iMpc^^ which already 
reaches ~ 50% at fc ~ 5/iMpc~^. The reason for these large 
differences is unclear. 

iGuillet et all (|2010l ) used the MareNostrum simulation, 
which w as run using the adaptive mesh refinement code 
RAMSES l|Tevssieij|2002l ') , to investigate the effects of baryons 
on both the variance and the skewness of the mass distri- 
bution. They used a 50 Mpc box with 1024^ dark mat- 
ter particles and included metal-dependent gas cooling, UV 
heating, star format ion, SN fee dback (using the kinetic feed- 
back prescription of [Dubois fc Tcvssicr 2008) and metal en- 
richment. AGN feedback was not considered. Unfortunately, 
they were not able to run their simulations down to z = 0, 
but quote results at redshift 2 instead. In order to better 
compare to their results, we have examined the power spec- 
tra of REFJL100N512 and DMONLY.L100N512 at z = 2. 
In our reference simulation the scale on which baryons sig- 
nificantly reduce the power increases with time (note that 
AGN shows the opposite behaviour, see Fig. [8]): in REF 
the 1% level is first reached at fc « 2/iMpc~^ for z = 2 
and at fe « 0.8/iMpc~^ for z = 0. Meanwhile, the effect 
on the power on scales fc > 10/iMpc~^ hardly changes, 
and the transition scale from a decrease to an increase in 
power relativ e to DM ONLY remains fixed at fc « 7 /i Mpc~^ . 
iGuillet et aLl(|2010l l. on the other hand, do not detect a sys- 
tematic decrease in power due to baryons at any scale. They 
find that the power is increased by 1% relative to a dark 
matter only simulation at fc ~ 3/iMpc~^, reaching 40% at 
fc « 10/iMpc-^ For our reference model we instead find a 
2% decrease for fc ~ 3/iMpc~^ and only a 6% increase at 
fc « 10 /i Mpc~^ . It is hard to say why these results lie so far 
apart, and especially why the baryons in their simulation do 
not reduce the power on large scales due to pressure effects. 

W e also compare to the recent study bv iCasarini et al. 
I2OI1I), who use the SPH code gasoline jWadslev et al. 
I2OOJ) to perform their simulations. They use two different 
volumes: a box of 64 Mpc on a side, and a much larger 
256 Mpc box, both with only 256'^ dark matter and an 
equal number of gas particles. Note that the mass resolution 
of their L = 64 /i"^ Mpc run is comparable to that of our 
fiducial run, while the resolution of their L = 256 Mpc 
run is much poorer. They include radiative cooling, a UV 
background, star formation and SN feedba c k. For the latter 
they use the prescription of IStinson et al.l (|2006l ). in which 
Type II SNe are modelled using an analytical treatment 
of blastwaves combined with manually turning off radiative 
cooling. Metal-line cooling and AGN f eedback were not con - 
sidered. Using their 64/i~^Mpc box. ICasarini et al.l (|201lD 
find an ~ 1% decrease in power at fc « 1 — 2/iMpc~^ 
and an increase in power at smaller scales, which reaches 
20% at fc 10/iMpc"^ These results are in reasonable 
agreement with both ljing et al . (2006) and our model NOZ- 
CO OL. However, when u sing their 256 Mpc box, they - 
like IGuillet et all (|2010l ) - find no decrease in power due 
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to baryons at any scale, but instead a steady increase in 
power that reaches 1% at A; ~ 1 — 3/iMpc~^ and 40% at 
fc«10/iMpc-\ 

Finally, we discuss the work by iLevine &: GnedinI 

l|2006 l). who used a toy model, rather than a hydrodynamic 
simulation, to evaluate the potential effect of AGN feedback 
on the matter power spectrum. In their models only the evo- 
lution of dark matter was followed explicitly. The gas was 
assumed to trace the dark matter at all scales and galaxy 
formation and the associated physical processes were not 
included. Their standard simulation volume is 64/i~^Mpc 
on a side, and the simulation was run with resolutions of 
1, 0.5 and 0.25 Mpc. We note that even their highest 
resolution is more than two orders of magnitude below the 
spatial resolution in our standard simulations. The gas was 
assumed to have a constant temperature of 1.5 x 10* K at all 
redshifts. A quasar luminosity function was used to deter- 
mine the number of AGN at a given redshift and luminosity, 
which were then each placed at a random location, although 
biased towards high-density regions. Of the AGN's bolomet- 
ric luminosity, a fraction ek = 1% was used to drive spheri- 
cally symmetric outflows. Within these outflow regions the 
baryon fraction was assumed to be zero. After computing 
the power spectrum, they found a large discrepancy between 
simulations with different resolutions: when using a resolu- 
tion of l/i~^Mpc, they found a reduction of roughly 10% 
in power for 0.3 < < S/iMpc"^ at z = 0, relative to a 
simulation which did not include AGN, while their higher- 
resolution runs produced instead an increase in power at all 
scales, of up to 20%. We found a decrease in power of 1% at 
k « 0.3/iMpc-\ reaching > 10% for 2 < fc < 50/iMpc-\ 
which does not agree with their results, even in terms of the 
si gn of the effect. Neverthe less, we do confirm the conclusion 
of Levine fc GnedinI (|2006l ) that AGN feedback can greatly 
affect the matter power spectrum on a wide range of scales. 

Even though our current understanding of galaxy for- 
mation still allows for significant deviations between studies, 
some qualitative results are the same: in the absence of AGN 
feedback, baryons will affect the matter power spectrum sig- 
nificantly on scales fc ~ 1 — lO/iMpc"^. Furthermore, all 
studies agree that the increase in power due to baryon s is of 
the order of 10% at fc = 10feMpc~\ Jing et all (120061'). ou r 



reference and NOZCOOL models, and lCasarini et al.l \2Qli ) 
for their high-resolution simulation all predict a relative de- 
crease in power of ~ 1% at fc ~ l/iMpc~^. iRudd et al.l 
(|2008l l also find a decrease in power due to baryons, but in 
their case the effect is far stronger than that of any other 
study, and is seen at much larger scales (fc < 1 /iMpc~^). 

However, like our reference simulation, all these simu- 
lations suffer from t he well-known overcool ing problem. As 
was demonstrated bv lMcCarthv et all (12010 *). the AGN sim- 
ulation does not. We have shown that the inclusion of AGN 
has a tremendous effect on the matter power spectrum for 
\ < 10 Mpc, both when compared to a simulation that 
includes only dark matter and when compared to simula- 
tions that include baryons and galaxy formation but not 
AGN feedback. Th erefore, contrary to what, for example, 
iGuillet et all (|2010h claim, simulations that suffer from over- 
cooling cannot be considered extreme models for which the 
effects of baryons on the total matter power spectrum are 
maximised. Instead, they are prone to tinrferestimate the ef- 
fects on large scales. Indeed, model AGN predicts a relative 



decrease in power of ~ 1% already at fc = 0.4/iMpc"^ The 
decrease in power reaches several tens of percent on scales 
fc ~ 1 — 10 /i Mpc~^, while simulations that suffer from over- 
cooling instead predict a strong increase in at least part of 
this range. Based on our results and on the comparison to 
other studies, we argue that the inclusion of AGN in cosmo- 
logical simulations is at present even more important than 
the improvement or convergence of existing prescriptions for 

other baryonic effects. 

Motivated by the results of iRudd et al.l (|2008h . 
IZentner et al.l (|2008l ) have proposed a method to account 
for the effects of galaxy formation on the matter power spec- 
trum. This method assumes that the effects of baryons can 
be captured by a change in the halo concentration-mass rela- 
tion. However, it is unlikely that such an approach can truly 
model the effects of baryons on the power spectrum. Since 
AGN-driven outffows significantly affect scales much larger 
than the sizes of individual haloes, the assumption made by 
IZentner et al.l ([2008) will certainly not be valid when AGN 
feedback is included. 



5 CONCLUSIONS 

Upcoming weak lensing surveys, such as LSST, EUCLID, 
and WFIRST aim to measure the matter power spectrum 
with unprecedented accuracy. In order to fully exploit these 
observations, theoretical models are needed that can pre- 
dict the non-linear matter power spectrum at the level of 
1% or better on scales corresponding to 0.1 Mpc < fc < 
10/iMpc~^. Here, we have employed a large suite of simula- 
tions from the OWLS project, as well as the highly accurate 
power spectrum estimator powmes, to investigate the effects 
of various baryonic processes on the matter power spectrum. 
These tools have also enabled us to examine the distribution 
of power over different mass components, the back-reaction 
of baryons on the CDM, and the evolution of the dominant 
effects on the matter power spectrum. 

Our most important finding is that the feedback pro- 
cesses that are required to solve the overcooling problem (i.e. 
the overproduction of stars), have a dramatic effect on the 
matter power spectrum. Such efficient feedback, most likely 
in the form of outflows driven by AGN, were not present 
in the simulations used in previous studi es of the effects 
of baryons on the matter power spectrum djing et al.l 20061 : 
iRudd et al.ll2008l : iGuillet et al1l2010l : ICasarini et all 1201 ll ). 
Although it was generally assumed that overcooling would 
make the simulations conservative, in the sense that they 
would overestimate the baryonic effects, we demonstrated 
that the opposite is true. The efficient outfiows that are re- 
quired to reproduce optical and X-ray observations of groups 
of galaxies, redistribute the gas on large scales, thereby re- 
ducing the total power by > 10% on scales fc > 1 hMpc~^ . 

We emphasise that the model from which we draw this 
conclusion, the simulation that includes AGN feedback, is 
not extreme. On the contrary, we consider it our most re- 
alistic model. iMcCarthv et al.l (|2010l . I2OIII ) showed that it 
provides excellent agreement with both optical and X-ray 
observables of groups of galaxies at redshift zero. In partic- 
ular, it reproduces the temperature, entropy, and metallicity 
profiles of the gas, as well as the stellar masses, star forma- 
tion rates, and age distributions of the central galaxies, and 



© 2011 RAS, MNRAS OOO.nifTSi 



Galaxy formation and the matter power spectrum 15 



the relations between X-ray luminosity and both tempera- 
ture and mass. 

We showed that metal-line cooling, star formation, and 
feedback from SNe all modify the matter power spectrum 
by > 1% on the scales relevant for upcoming surveys. In 
the absence of AGN feedback, the simulations with baryons 
have ~ 1% less power relative to a dark matter only simu- 
lation on scales 0.8 < A: < 6/iMpc~^ (a consequence of gas 
pressure) and > 10% more power for k > lO/iMpc"^ (a con- 
sequence of gas cooling). However, as we noted above, AGN 
feedback can decrease the power for ^ k < 10/iMpc~^ 
by up to several tens of percent. Furthermore, some imple- 
mentations of stellar feedback, e.g. the strong SN feedback 
resulting from a top-heavy stellar initial mass function in 
starbursts, can create differences of the same scope and mag- 
nitude by redistributing gas out to very large scales. The 
effects from such baryonic processes on the matter power 
spectrum can even exceed those of a very large change in 
cosmology (e.g. WMAP3 to WMAPl). Indeed, differences 
> 1% persist even up to scales as large as those correspond- 
ing to fc « 0.3/iMpc"\ 

In the absence of AGN feedback, the back-reaction of 
baryons on the dark matter increases the power in the CDM 
component by 1% at ~ 2/iMpc~^ and the effect becomes 
larger towards smaller scales. However, when AGN are in- 
cluded they redistribute sufficiently large quantities of gas 
out to large radii to lower the power in the dark matter com- 
ponent b y 1 - 10% for 3 < fc < 30/iMpc"\ This is consis- 
tent with iDuffvet al.l(|2010l ). who used the same simulation 
to show that AGN feedback decreases the concentrations of 
dark matter haloes of groups of galaxies. We stress, how- 
ever, that the back-reaction of AGN feedback on the CDM 
will not be straightforward to implement in dark matter only 
models. While it may be possible to roughly model the effect 
of baryons in simulations without efficient feedback by rais- 
ing t he concentrat ion parameters of the dark matter haloes 
fe.g. IZentner et al . 2008), feedback from AGN redistributes 
the gas on scales that exceed those of their host haloes. 

The difference between dark matter only simulations 
and simulations that do include baryons is nearly the same 
for the WMAP3 and WMAP7 cosmologies, at least when 
AGN are included. This suggests that the relative effect of 
the baryons is roughly independent of cosmology, which will 
simplify future studies aiming to disentangle the two. 

For our most realistic simulation, which assumes the 
WMAP7 cosmology and includes AGN feedback, the differ- 
ence in power relative to the corresponding dark matter only 
simulation does not evolve much for z < 1 on large scales 
ffc < lOfeMpc"^). This is consistent with iMcCarthv et al.) 
1I2OIII ), who showed that the expulsion of gas through AGN 
feedback occurs mostly at z ~ 2 — 4, in the progenitors of 
today's groups and clusters of galaxies. 

We demonstrated that our conclusions are robust with 
respect to changes in the size of the simulation box and 
changes in the resolution (see Appendix A), with any addi- 
tional modelling uncertainties only making it less likely that 
the matter power spectrum can be predicted with 1% accu- 
racy any time soon. Looking at the large differences that still 
exist between the results of different authors, it is clear that 
much work remains to be done in understanding processes 
such as gas cooling and outflows. 

In a forthcoming paper (Semboloni et al., in prepa- 



ration), we will study the implications of our findings for 
weak lensing surveys in more detail. In this work we will 
also demonstrate that the use of optical and X-ray observa- 
tions of groups of galaxies can significantly reduce the un- 
certainties in the predictions of the matter power spectrum. 
While this provides a strong incentive for obtaining better 
and more observations of groups of galaxies, it is important 
to note that such auxiliary data will never completely re- 
move the uncertainty inherent to cosmological probes of the 
matter distribution on scales that are potentially affected 
by baryonic physics. This is because one can never be sure 
that all the relevant effects are constrained by the secondary 
observations. For example, it may be that other models for 
the baryonic physics exist that also reproduce optical and 
X-ray observations of groups, but nevertheless predict dif- 
ferent power spectra. It will therefore be crucial to consider 
a wide variety of observations, with optical and X-ray as 
well as Sunyaev-Zel'dovich observations holding particular 
promise, and a large range of models. 

While the strong baryonic effects that we find imply 
that the cosmological constraints provided by upcoming 
weak lensing surveys will be model-dependent, it also means 
that such surveys will provide constraints on the physics of 
galaxy formation on scales that are difficult to obtain by 
other means. 

Tabulated values of power spectra for redshifts z — — 
6 are available for all the simulations shown in this paper 
as Supporting Information with the online version of this 
article, and at http : //www. strw. leidenuniv .nl/VDll/ (see 
Appendix B). 
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APPENDIX A: CONVERGENCE TESTS 

Here we investigate the effects of changing the box size or 
resolution of the reference simulation on its power spectrum. 



Al Box size 

In Figure lAll we vary the size of the box at constant res- 
olution. The difference between the power spectrum of the 
100 and the 50/i"^Mpc box is smaller than the difference 
between the latter and the 25/i~^Mpc box for all k, and 
the power spectrum of the largest box is nearly converged 
for k > 20/iMpc-^ 

However, there are differences of up to a factor of a few 
at larger scales. For reference, we also show the input power 
spectrum linearly evolved to z = a nd the HALOFIT m odel 
of the non-linear power spectrum bv lSmith etahl (|2003l ) (see 
i]3.1|) . The first wave mode corresponds to the size of the sim- 
ulation box, which means that the power measured on this 
scale is meaningless; hence, we have omitted this point in 
all of our figures. The second and third wave modes closely 
follow the linear power spectrum. Note that the curves have 
very similar shapes on large scales, with the larger boxes 
shifted to larger scales. This is a consequence of employing 
the same seed for the random number generator used to cre- 
ate the initial conditions. Perturbations that should go non- 
linear (A < 10 Mpc) are unable to collapse if their wave- 
length is close to the size of the box, which in turn suppresses 
the power on smaller scales. One might therefore worry that 
even the 100 Mpc box is not large enough to obtain ac- 
curate power spectra for k > I ^Mpc~^. Lacking larger sim- 
ulations to check this, we compare to the HALOFIT model 
for the non-linear power spectrum, which shows where the 
transition from the linear power spectrum should take place 
at redshift zero. The power spectrum for the 100 Mpc 
box follows this model very well on large scales, suggesting 
that a simulation of this size is very close to converged. 

Note that finite volume effects only prevent us from ob- 
taining highly accurate absolute power spectra, and only for 
the largest scales, while our results are based on the rela- 
tive comparisons between models that used identical initial 
conditions. Since the 100ft~^Mpc box extends up to the 
largest non-linear scales, and since all simulations start from 
the exact same realisation of the linear power spectrum at 
z = 127, we do not expect our results to be affected by the 
finite volume of the simulations. 



A2 Numerical resolution 

In Figure IX2I we investigate the effects of changing the reso- 
lution for the reference simulation by varying the number of 
particles while keeping the box size fixed. The power spec- 
trum of REF.L100N128 is quite noisy for k > lOO/iMpc"^ 
because of its much higher Poisson noise level. Testing 
for convergence on these scales is only possible thanks 
to the accurate shot noise subtraction. Surprisingly, the 
relative difference in power between REF_L100N128 and 
REF_L100N512 is smaller than the difference between the 
latter and REF_L100N256. When increasing the resolution 
beyond that of REF_L100N256, the power begins to de- 
crease. To examine if this trend continues, we compare the 



\ [Mpc/h] 

100.00 10.00 1.00 0.10 0.01 




0.1 1.0 10.0 100.0 1000.0 

k [h/Mpc] 

Figure Al. Test of convergence of the 2 = matter power spec- 
trum in the reference model with respect to the size of the simu- 
lated volume, where the box size and particle number are varied 
in such a way as to keep the resolution constant. Also shown are 
the linear input po wer spectrum and the analytical non-linear 
power spectrum by ISmith et al.l l|2003h . The red, dotted line in 
the top panel shows the (subtracted) theoretical shot noise level. 
The bottom panel shows the ratio of REF_L100N512 with respect 
to the other simulations. 



power spectrum of REF_L050N256, which has the same res- 
olution as REF.L100N512, to that of REF.L050N512 in the 
panel on the right. We see that the power on the smallest 
scales (k > 10/iMpc~^) converges only slowly, but that the 
trend of decreasing power with increasing resolution contin- 
ues. This may indicate that, as lower mass haloes become 
resolved, the overall effects of supernova feedback become 
stronger. 

We can verify this by isolating the effects due to baryon 
physics from those due to a more straightforward depen- 
dence on resolution. To this end, we examine what the effect 
is of increasing the particle number of the DMONLY simu- 
lation with a 100 Mpc box in Figure fXSl The behaviour 
here is quite different: as A'^ grows, more low-mass haloes 
are resolved and the power on small scales increases. As we 
observe a reversed trend in Figure [X2I we conclude that the 
increased baryonic effects that accompany a higher particle 
number are more important for the power spectrum than 
the straightforward dependence on resolution. 

The difference between REF_L050N256 and 
REF.L050N512 is ~ 0.1% at fc = 1/iMpc"^ and 
~ 2% at fc = lO/iMpc"'^. We conclude that simulation 
REF_L100N512 is sufficiently converged for the scales of 
interest for this study, k < 10/iMpc~^. Note that, since 
we are only interested in the relative differences between 
simulations with equal resolution, the uncertainty will in 
practice be much smaller. With increased resolution we ex- 
pect feedback processes to become more effective, meaning 
that we may have underpredicted the differences between 
models with different feedback processes in low-mass haloes 
on small scales. 

Similar tests were performed by IColombi et"al] (|2009l ') 
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Figure A2. As in Figure lAl] but now the numerical resolution is varied while keeping the box size constant. The left and right panels 
show power spectra for box sizes of 100 and 50 Mpc, respectively. On small scales, k > lO/iMpc"^, convergence is slow. 
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Figure A3. Same as the left-hand panel of Figure lX2l but now for 
DMONLY instead of REF. Here the behaviour is as expected: as 
the number of particles goes up, more low-mass haloes form and 
the power on small scales increases. A comparison with Figure[A2] 
shows that increasing the resolution leads to stronger baryonic 
effects which may reverse the sign of the trend with resolution. 



for the convergence of powmes, which keeps the statistical 
error bounded through its use of foldings. Its value depends 
on the quantity C{k), which is defined as the number of 
independent wave modes at a given wave number k; to be 
more precise, we ap proximately have AP/P oc C{k)~^^^ 
IColombi et al.l 120091 ). For our fiducial grid with 256^ grid 
cells, one can expect the statistical error to remain below 
|AP|/P ~ 1.2% as long as errors due to shot noise do not 
dominate. We have checked that this is indeed the case. Note 
that this means that we can confidently measure 1% diflFer- 
ences between simulations using our fiducial values, as we 



Table Bl. Power spectrum values for AGN.WMAP7.L100N512 
for a subset of scales at z = (full table available online). 



z 


k [h/Mpc] 


P(fc) [h-^Mpc^] 


A2(fc) 


0.000 


0.12566371 


4364.4776 


0.43876514 


0.000 


0.18849556 


1853.4484 


0.62886024 


0.000 


0.25132741 


1524.3814 


1.2259802 


0.000 


0.31415927 


1112.5603 


1.7476056 


0.000 


0.37699112 


847.62970 


2.3007519 



are interested in systematic offsets covering at least a small 
range of scales in fc-space, rather than random deviations. 



APPENDIX B: TABULATED POWER 
SPECTRA 

Table IBll shows the power spectrum values for 
our most current and realistic simulation to date, 
AGN_WMAP7.L100N512, for a subset of scales at z = 0. 
Our fiducial powmes values of 256^ grid points and 7 fold- 
ings were used, and shot noise has been subtracted. The full 
table, with power spectrum values at all scales shown in this 
paper and redshifts up to 2: = 6, as well as tabulated data for 
all other simulations presented in this paper, are available 
as Supporting Information with the online version of this 
article and at http://www.strw.leidenuniv.nl/VDll/. 

This paper has been typeset from a Tp^X/ ET^rjX file prepared 
by the author. 
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